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ABSTRACT 

This paper is concerned with the theoretical properties of high contrast coronagraphic images 
in the context of exoplanet searches. We derive and analyze the statistical properties of the 
residual starlight in coronagraphic images, and describe the effect of a coronagraph on the speckle 
and photon noise. Current observations with coronagraphic instruments have shown that the 
main limitations to high contrast imaging are due to residual quasi-static speckles. We tackle 
this problem in this paper, and propose a generalization of our statistical model to include the 
description of static, quasi-static and fast residual atmospheric speckles. The results provide 
insight into the effects on the dynamic range of wavefront control, coronagraphy, active speckle 
reduction, and differential speckle calibration. The study is focused on ground-based imaging 
with extreme adaptive optics, but the approach is general enough to be applicable to space, with 
different parameters. 

Subject headings: instrumentation: adaptive optics, instrumentation: high angular resolution 



1. Introduction 

Direct imaging of faint companions or planets around a bright star is a very difficult task, where the 
contrast ratio and the angular separation are the observable parameters. The problem consists of detecting 
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a faint source over a bright and noisy background, mainly due to the diffracted stellar light. High contrast 
ratios and small angular separations correspond to the most difficult case. Typically, for extra-solar giant 
planets, contrast ratios of about 10~ 7 are exp ected in the near infrared (J,H,K bands), based on models 



for re latively young objects of about f 00 Myr (jChabrier fc BaraffduOOOt iBaraffe et al.ll2003t iBurrows et al 



20041 ). According to these models, older objects would be an order of magnitude fainter. Terrestrial planets 
are much fainter than giant planets, about 3 to 4 orders magnitudes fainter depending on the wavelength 
range. 

In the case of ground-based observations with Adaptive Optics (AO), the residual uncorrected aberra- 
tions produce random intensity fluctuations of the background, which appear as speckles in the field. In direct 
non-coronagraphic high quality images, these speckles mainly appear at the position of the diffraction rings 
of the sta r. In Figffl we show two high-quality po int spread functions (PSF) obtained with the AO system at 
Palomar (jTrov et al.ll2000t lHavward et al.l 120011). where the speckles are clearly visible. This phenomenon, 
also known as "speckle pinning" feloemhof et ak 2001), can be explained us i ng an expansion o f the point 



spread function (PSF) (jsivaramakrishnan et al.ll2002l . |20ol lBloemhoj|2003at [perrin et alJliool iBloemhoi 



2004b ). An alternative approach using a statistical model ( Aime fc Soummer 2004bh enables a deep 



derstanding of the phenomenon, especially from the statistical point of view, by providing information on 
the speckle variance. The effect of a coronagraph on speckle noise is well explained with this approach. In 
particular, we show that static aberrations produce residual speckle pinning after a coronagraph, which has 
important impl ications in high contra st imaging. Statistical properties of long exposure AO images were 
also studied by Fusco fc Conanl ( 2004 ). 



The dynamic range of a coronagraphic image corresponds to the faintest companion that can be detected 
at a given position in the field, at the detection limit. It is usually expressed as a magnitude difference or an 
intensity ratio, relative to the unocculted central star at some signal to noise ratio (SNR) level. Although 
the dynamic range is a two-dimensional map of the detection sensitivity in the field, it is often represented 
as a radial profile showing the magnitude difference Am as a function of the angular separation r. A radial 
profile should be acceptable for most applications, especia lly when the coronagraph does not have particular 
asymmetries (|Lvotill939tlAime et al.ll2002l ; ISoumm cr 2005). Otherwi se, a two-dimensional map of sensitivity 



is required when the co ronagraph shows an asymmetric response (|Kuchner fc Spergell l2003t iKasdin et al 



2003c iRouan et alJ 120001 ) . Understanding, measuring and predicting the dynamic range is still one of the 



important issues in this field, with implications for instrumentation (design, observing strategies) and data 
reduction and analysis. 

During the design phase of future planet finder instru ments using Extreme Adapt ive Optics (ExAO) 
and coron agraphy, for examp l e the Gemini Plan et Imager ([Macintosh et al.ll2004i bood ). or the ESO/VLT 



SPHERE (Beuzit et al 



2005 



Fusco et al.ll2005l ). it is necessary to anticipate what part of the observable 



parameter-space (contrast vs. separation) can be probed, and link it to the actual p hysica l parameter space 
(mass vs. semi major axis). Such studies have been carried out by Brown ( 2004bM a. 2005h in the context of 
terrestrial planet searches. 



When operating an existing high-contrast instrument, like the Lyot pr oject coronagra 



2004 ) , the dynamic range has to be measured to evaluate the perfor mance ( Soumme r et al 



ph (IQ ppcnhci mer et al 



2006; Hinkley et al 



2006) . In the case of a dete cted object, photometry and astrometry ( Digbv et al . 2006; ISivaramakrishnan fc Qppenheimer] 



20061: iMarois et alj|2006bh are necessary to help determine the objects characteristics. Dynamic range com- 
putations are also important in the case of non-detection, to determine which part of the parameter space 
has been probed by the experiment, and which physical objects can be ruled out. 
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The dynamic range is limited by the intensity fluctuations close to the star. These are due to sev- 
eral sources: speckle noise, photon noise, detector noise, background noise etc. Speckle no i se is known 



to be the main limitation for high cont rast imaging, either in di rect images (Marois et al 



Masciadri et al 



2001 



2005; Marois et al 



Boccaletti et al 



2003 



2006 a|), or coronag raphic images (jBeuzit et al.lll997l ; lOppenheimer et al 



2003 



2005 



2004t iHinklev et al.ll2006l ) . Speckles find their origin in wavefront imperfections 



(amplitude and phase errors), whether they correspond to uncorrected atmospheric residual errors (residual 
atmospheric speckles) or slow-varying wavefront caused by mechanical or thermal deformations. The main 
pro blem comes from th e se quasi-static speckl e s, which can be calibra ted either using active pr e-focal meth- 
ods 



Malbet et al 



2000: 



( 19951) :lGive'on et all (|2006l);|Borde fc Traubl ( 2006 ) or using post-processing (jMarois et al 



Sparks fc Ford 



has been achieved by 



2002; Marois et al. 2006a|). An alternative approach based on non- redundant masking 
Llovd et al. ( 20061 ). Assuming a good enough calibration of these quasi-static speckles, 
the physical limit of the system is set by the residual atmospheric aberrations. The noise limitation in 



high dynamic range image s has been studied by several authors, using n umerical simulations ( Boccaletti 



2004; 



Cavarroc et al 



Perrin et a l. 2003; Guvon 



2005; Sivaramakrishnan et al 



20051) 



2006h. or other theoretical or empirica l approaches (jAngell 1994 iRacine et al.lll99*9t 



In lAime fc Soummerl (|2004bl ). we modeled the statistics of AO-corrected, direct images, and discussed 
the effect in terms of signal to noise ratio on these images qualitatively. In this paper, we examine the effects 
a coronagraph has on the properties of residual speckles. A semi-analytical method to compute the dynamic 
range can be derived from the statistical properties of the speckle and photon noise. We compare our re- 
sults with purely numerical simulatio ns. Results in th i s paper are presented using an Apodiz ed Pupil Lyot 



Coronagraph (APLC) as an example (lAime et al 



2001 



200 



2 



Soummer et al 



2003a 



Soummer 



2005), but are valid 



for any other type of mask coronagraph ( Rouan et al. 2000l: iKuchner fc Traubll2002l: [Soummer et al.ll2003_ 
Mawet et al .l2005h. pure apodizers or shaped pup i ls Jjacquinot fc Roizen-Dossierll964lNisenson fc Papaliolio; 



Soummer et al. 2003ri Kasdin et al. 20031 Aime 2005b ). Furthermore, the statistical model can be 



modified to include both static and quasi-static aberrations, and we discuss the coherent interaction between 
residual atmospheric and quasi-static aberrations. 

Our theoretical model and results apply for both space and ground-based imaging. However, we illustrate 
the results with simulations in the case of ground-based Extreme Adaptive Optics (ExAO) and coronagraphy. 




Fig. 1. — Direct images of the bright star HD 137704 (magnitude V=5.47) obtained with the Adaptive 
Optics system at the Palomar Hale telescope on June 6th 2004. The core of the image has been saturated 
to illustrate the speckles. The bright speckles localized at the position of the diffraction rings are called 
pinned speckles. Their statistical properties in direct and coronagraphic images are discussed at length in 
this paper. 
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2. Propagation through a coronagraph in the presence of aberrations 



In this section we assume that all the static aberrations, mainly optical quality (polishing for example) 
and misalignment errors, c an be represen t ed at t he entrance pupil of the telescope. The case of out of pupil 
aberrations is discussed by iMarois et al.l (|2006cf ) and does not affect our monochromatic statistical model. 
We consider a n instrument with an ExAO system and a generic coronagraph which can describe any type of 



mask designs (IRouan et al 



2000 



Aime et al 



2002; 



,nd a generic co ronagrap 
Kuchner fc Traubf 2002: 



Soummer et al 



2003c; 



Soummer 



2005t iMawet et al.l 120051 ). The formalism als o applies to the sh aped pupil approach which corresponds to 
the case of direct imaging with apodization ( Kasdin et al. 20031) . The ExAO and telescope characteristics 
are chosen to be consistent with current or future projects on eight-meter telescopes and illustrations are 
given for an APLC. 

These coronagraphs consist of optical filtering in four successive planes denoted 1,2,3,4 hereafter. The 
first plane corresponds to the entrance aperture (possibly apodized), the second plane is the focal plane 
where a mask is applied (opaque, graded or phase-shifting), the third plane corresponds to a relay pupil 
plane where a diaphragm is applied (the Lyot Stop), and finally the fourth plane corresponds to the final 
focal plane. 

In the pupil plane, we model the wavefront using a static phase aberration </? s (r), a residual random 
atmospheric phase ^(r), and amplitude aberrations /c(r), including the eventual apodization. The complex 
amplitude is: 

*i(r) = P(r) p(r)e j ^ + ^\ (1) 
where the function P(r) describes the normalized aperture transmission: J P(r)dr = 1, and r = (x, y) i s the 



coordinate vector, used in both pupil and field. Following the notations of I Aime fc Soummer ( 2004b ). we 
can write the wavefront complex amplitude at the entrance pupil as the coherent sum of three terms: 



*i(r) - [A + A s (r) + a(r)] P(r), 



(2) 



where A is a deterministic term corresponding to a perfect plane wave, ^4 s ( r ) is a deterministic complex 
term corresponding to the static aberrations, and a(r) is a random term with zero mean (E[a(r)] = 0) 
corresponding to the uncorrected part of the wavefront. The probability density function (PDF) of this 
complex amplitude is illustrated in Fig[2] without and with AO. Static aberrations are not included in this 
figure. A is defined as the mean of the complex amplitude, averaged over the pupil: 



A = E[J *i(r)P(r)dr]. 



(3) 



\A\ 2 is therefore the Strehl Ratio of the system (Hardy 19981) . 

With E[a(r)] = 0, Eq[2]implies E[Wi(r)] = A + A s (r). Integrating this equation, we obtain: 



E[f Vi(r)P(r)dr]=A + J A s {r)P{r)dr, 



and therefore: 



J A s (r)P(r)dr = 0. 
Assuming that the phase errors are stationary over the aperture, we obtain: 



a(r)P(r) » 0. 



(4) 
(5) 

(6) 
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Fig. 2. — Probability Density Functions (PDF) of the pupil plane complex amplitude. Left: PDF of un- 
corrected atmospheric wavefronts obtained from a von Karman power spectrum model, with an outer scale 
Lq = 20m and a seeing loq — 0.8 arcsec. The phase excursion is uniform over (0 — 2tt) and the thickness of 
the annulus corresponds to the amplitude scintillation. Center and Right: PDFs of a AO-corrected wave- 
fronts for Strehl Ratio of 90% (center) and 95% (right). These distributions in the complex plane look like 
de-centered crescents. The length of the crescent corresponds to the phase excursion, and the thickness to 
scintillation, which is assumed to have no effect on the AO system. The effect of the improved AO correction 
on the phase excursion is obvious between these two figures. 



Note that the difference between Eqj5]and Eqj6] comes from the fact that A s (r) is deterministic and a(r) is 
random: Although A s (r) can be defined with zero-mean over the aperture, each independent realizations of 
a(r) do not necessarily have an exactly zero average over the aperture. 

The specific case of quasi-static aberrations is not considered in this section and will be treated in Sec l3.2l 
The three terms of Eq|2] are illustrated in Fig|3l The length of the vector A is arbitrary in the figure, to 
illustrate that the modulus of A is not unity and that the vectors A, A Sl and a are defined according to the 
definitions above. In the first focal plane, a coronagraphic mask is applied at the center of the image of the 




Fig. 3. — Illustration of the decomposition of the wavefront as the sum of three complex vectors. We consider 
a static phase term ip s and a residual atmospheric phase ip. A is the mean wavefront and \A\ 2 is therefore 
the Strehl Ratio. A s corresponds to the static aberrations, and a to the zero-mean error term. 

star. Writing the mask transmission as 1 — M(r), allows us to accommodate any type of mask coronagraph, 
including Lyot, APLC, Band-Limited, Phase Masks. For example, a classical hard-edged Lyot coronagraph 
(or APLC), is described using a top- hat function for M. The complex amplitude of the wave in the focal 
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plane is given by a scaled Fourier Transform (FT) of this pupil amplitude (|Goodmanlll996l ): 

* 2 (r)=^[*i(r)](l-M(r)), 



(7) 



where the symbol T denotes the scaled FT. For clarity, we will omit the wav e length-dependen t scalin g 
factors. For t he complete chromatic formalism see for example lAime et al.l (|2002l) ; ISoummer et al.l (|2003a ) ; 
Aimel ( 2005a ). In the next pupil plane, the complex amplitude before the Lyot stop P'(r) is also the sum of 
three terms: 

* 3 (r) = A* c (r) + * s (r) + * Q (r), (8) 



w here ^(r) 



P(t) — P (r) * T[M (r)] is the complex amplitude in the lyot stop plane for a perfect wavefront 



(jSoummer et al.ll2003al ). The two other terms ^(r) and ^ a ( r ) correspond to the propagation of the terms 
A s (r) and a(r) respectively. For example: 



*„(r) = a(r)P(r) - a(r)P(r) * T[M{v)\. 



(9) 



The perfect coronagraph term \I/ C (r) and the term \t a ( r ) are shown in intensity in FigUJ The coronagraphs 
rejects most of the starlight outside the image of the aperture in the Lyot plane for the perfect part of the 
wave, but most of the energy remains inside the aperture for the speckle part. The Lyot stop is applied 




Fig. 4. — Example of propagation though a coronagraph for a monochromatic APLC, corresponding to Eqj8] 
Left: intensity in the Lyot stop plane in the perfect case (|^ c (r)| 2 ) . Most of the light is rejected outside of 
the pupil aperture and will be eliminated by the Lyot Stop. Right: intensity of the speckle term |\E' Q (r)| 2 , 
assumed alone. Most of the energy remains inside the geometric aperture and will not be eliminated by the 
Lyot Stop and appear as residual speckles in the final image. Both images are represented with the same 
scale. 

in this plane. In the case of an APLC, the Lyot stop is identical to the entrance pupil; in all other cases, 
the Lyot stop is undersized. With P(r)P'(r) — P'(r), and with the notations S(r) — !F[a(r) P'(t)] and 
S s (r) = T[A s {y) P'(r)], we obtain the complex amplitude in the final focal plane: 



* 4 (r) = AV d (r) 

+S s (r) - (S s (r)M(r)) * T[P'(r)} 
+S(r)-(S(r)M(r))*.F[P'(r)], 



(10) 



where ^>d d e notes the focal wave amplitude of the coronagraph in the perfect case, following the notations of 
Aime et al.l l|2002h . The convolution product (S(r)M (r))*.F[P'(r)] in EqIT0]has a negligible effect outside the 
mask area. Indeed, the spatial extension of (S(r)M(r)) is limited to the occulting mask area, and J-[P'(r)] 
is a rapidly decreasing function (the Airy amplitude in a perfect case) whose characteristic size is X/D. The 
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result of the convolution does not extend much beyond the mask area, which is illustrated in FigJSl where 
we show an example of the effect of the convolution term using a numerical simulation. This can also be 
explained by considering the propagati on of phase ripples through a coronagraph and constructing a Bode 
diagram (jSivaramakrishnan et al.ll2007l ). 



Fig. 5. — Illustration of the terms of EqfTOl Left: modulus of the speckle term alone S(r) = T[a{v) P(r)]. 
Right: actual term "^(r) represented in modulus at the same scale. The cross terms of Eq[lO]only affect 
the central part of the image, inside an area corresponding to the occulting mask size. 



Grouping all the deterministic terms together, we introduce 

C*(r) = A * d (r) + S s (r) - (S.(r)M(r)) * F[P'(r)], 



(11) 



and we obtain a similar expression to the case without coronagraph where C(r) = A jF[P(r)] (jAime fc Soummer 
2004bh : 

* 4 (r) = C(r) + S(r). (12) 

The wave amplitude after a coronagraph appears as a sum of a deterministic term C(r), and a random term 
S(r), at each position in the focal plane, outside the mask area. The deterministic term C(r) corresponds to 
the focal plane complex amplitude of the coronagraph in the presence of static abe rrations. The random term 
S(r), associated with the speckles, is identical to the case without coronagraph (jAime fc Soummerll2004bl ): 
the coronagraph has a negligible effect on the random part of the wavefront, as illustrated in Figf4]and Figl5] 
Formally, the effect of the coronagraph is to replace the wave amplitude without coronagraph C(r) by the 
coronagraphic amplitude C (r) . In the case of pure apodizers (shaped pupils) , the direct apodized term C (r) 
is used. 

The random term S(r) is non stationary in the field. The profile for S(r) can be computed from a 
simulation of the AO system, as we detail in Sec l4.21 Low-order a berrations can also be in c luded in thi s 
description, but usually require a specific study, as for example in ISivaramakrishnan et al.l (|2005l . 12007). 
The profile or the two-dimensional map for C(r) can be computed independently, considering a perfect 
coronagraph in the presence of deterministic static aberrations (and normalized to the SR). Even in the case 
of an ideal coronagraph that cancels all the star light for a perfect wave, the deterministic term C(r) still 
contains the terms due to the static aberrations and will contribute to speckle pinning, as discussed below. 
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3. Statistical properties of direct or coronagraphic images 
3.1. Statistical model and properties of speckles 

3.1.1. Complex amplitude distribution 

In this section, we discuss the distribution of the complex amplitude in the focal plane. We consider the 
case of monochromatic direct images for simplicity. The case of coronagraphic images is formally identical to 
the coronagraphic case, according to the approximations described in the previous paragraph fEq J12|). The 
focal plane complex amplitude is the Fourier transform of the pupil plane complex amplitude: 

* 2 (r) = J P(u) {A + a(u)) e~ 247rU r du. (13) 

The complex amplitude in the focal plane is therefore a sum of the random complex term a(r) weighted by 
the Fourier complex phasors. At the center of the image, the Fourier ph asors vanish, so a special treatment 
for this particular case (and the transition region around it) is necessary ( Soummer fc Ferrarll2007 ). Outside 



the central point of the image, the distribution of the complex amplitude can be derived using known results 
in Signal Processing, based on reasonable assumptions. We assume that the complex amplitude in the pupil 
plane can be represented by discrete values (an implicit assumption in any numerical simulation), and that 
the correlation of the complex amplitude between two points in the pupil plane decreases with distance 
between them. Under these hypotheses, it can be s hown that the d istribution of the complex amplitude in 



the focal plane is asymptotically circular Gaussian (jBrillingerl 1 1 98 lh . We remind the reader here that if the 



real and imaginary parts of a complex number z are Gaussian, its distribution is said to be Gaussian. If 
the real and imaginary parts are independent and have same variance, the d istribution is said to be circular 
Gaussian and denoted z ~ A/" c (0, a 2 ). See Fig. 2 of Aime fc Soummer ( 2004bl ) for an illustration of the focal 



plane PDFs. The circularity of the Gaussian distribution is due to the Fourier phasors mixing the complex 
amplitudes in the complex plane in the Fourier integral (|13[) . where u varies between —D/2 and D/2. For 
positions r in the focal plane such as r > X/D, the Fourier phase term therefore varies between and 2-7T and 
this circularization occurs. The complex amplitude of the wave in the focal plane ^(r) follows a circular 
Gaussian law, de-centered by the mean of the amplitude C(r) and denoted: ^(r) ~ Af c (C(r), E[|5(r)| 2 ]). 

In FigEJ we give an illustration of the distribution of the complex amplitude in the four successive planes 
of the coronagraph. This illustration is based on numerical simulations of a perfect APLC coronagraph and of 
an ExAO system. In the first pupil plane, we have a de-centered crescent (see Fig|2]). In the first focal plane, 
in this example at the top of an Airy ring, C(r) has a high absolute value, and the distribution is Gaussian, 
de-centered by this amount. Detail ed illustrations of the dece ntered Gaussian statistics as a function of the 
position in the field can be found in Aime fc Soummer ( 2004bh . In the following Lyot plane the coronagraph 



almost completely removes the perfect part of the wave (see FigSJ) , and the resulting distribution is similar 
to the initial distribution of the complex amplitude in the pupil, but centered at the origin. Finally in the last 
focal plane, without static aberrations, C(r) ~ 0, as — and the result is a centered circular Gaussian 
distribution. 



3.1.2. Intensity distribution 



In this section we derive the Probability Density Function (PDF) of the intensity from that of the wave 
complex amplitude. Our problem is formally equivalent to the case of laser speckles added to a coherent 
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Fig. 6. — Complex probability distributions in the four successive coronagraphic planes (pupil, focal, pupil, 
focal), at an arbitrary angular position in the field (r = 2.6X/D). This simulation corresponds to the case 
of an APLC without static aberrations. The distribution in the pupil plane corresponds to a typical ExAO 
for an 8-m class telescope delivering 90% SR, and including scintillation effects. In the focal plane, the 
distribution at a given position is a decentered Gaussian. In the Lyot stop plane after and APLC, the 
coronagraph has removed the deterministic part of the entrance pupil wavefront. In the final focal plane, 
the distribution is close to circular Gaussian distribution. 



background, which has been studied extensively (jGoodmanl Il975l 120061 ) . in particular in the context of 
holography. We introduce the two intensity terms: 



Ic = |C(r)| 2 
Is = E[|S(r)| 2 ] 



(14) 



Note that I c and I s are both functions of r, and that I c can describe both the direct or coronagraphic 
case, with and without static aberrations. Following Goodman, the joint PDF for the intensity and phase 
can be obtained from PDF of the complex amplitude, using the simple cartesian-polar change of variables 
(rj = \fl cos[0], £ = vTsin[#]), where the modulus of the Jacobian of this transformation is ||J|| = 1/2 and 
integrating the phase 8 to find the PDF for the intensity. 

An alternative derivation of the PDF for the intensity is to consider the properties of Gaussian dis- 
tributions. As discussed in the previous section, the speckle term S(r) is a circular Gaussian distribution 
S(r) ~ A/" c (0, I s ). The instantaneous intensity corresponding to the complex amplitude of Eqll2l is simply: 



/ = |S(r)+C(r)| 2 
= tar) + S(r)] 



($t[C(r) + S(r)\ 



(15) 



where 3? and 9 denote the real and imaginary parts. Using the properties of circular Gaussian distributions, 
9?[(7(r) + S(r)} and Q^C^r) + S(r)] are independent Gaussian random variables of same variance J s /2. We 
can rewrite the intensity with real and imaginary terms of variance unity: 



I = 



I, 



(jilv/^CXr) + 5(r)]) 2 + (^[V^Cir) + S(r)] 



2 ' 



(16) 



where var[K[\/27f T C , (r) + S(r)\] = var[S[v / 2if T (7(r) + S(r)]\ = 1. 

The random variable / follows a de-centered \ 2 with two degrees of freedom: xl( m ); w itfi a decentering 
parameter m = 2I^ 1 I C , ([Johnson et al.lll995l chap. 29). The probability density function for / is therefore: 



P(v) 



(17) 
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where f q (z) is the regularized confluent hypergeometric function and qFi(; q; z) the confluent hypergeometric 
function defined as: 



m = E 



i 



T(q + n)n\ 



oFi{;g;z) 

r(g) 



(18) 



Finally, the probability density function of the intensity I = I S /2I is: 



Mi) 



(19) 



This expression is equiva le nt [j] to the modified Ricia n distribution derived by 



Goodman 



( 1975) and used by 



Caeigal fc Canalesl (|l998h : ICanales fc Caeigall (|l999h : ICagigal fc Canalesl (EoO0l) : ICanales fe Caeigall d200lh : 



PI (I) 



exp 



(20) 



where Io denotes the zero-order modified Bessel function of the first kind. The Rician distribution is illus- 
trated in FigJT] A comparison between the Rician model and simulation data will be presented in Scc l3.1.4l 



An interesting particular case is when the background C(r) is zero and only the speckle term is present. 
Making I c = in Eq(20] (this happens at the zeros of the perfect PSF or using a perfect coronagraph), the 
PDF reduces to: 

M/) = iexp(-£). (21) 

This PDF corresponds to t he well-known ne gative exponential density for a fully developed speckle pattern 
(e.g. laser speckle pattern ) (lGoodmarJl2000h . Finally, the distribution at photon counting levels can be ob- 
tained performing a Poisson-M andel transformation of th e high flux PDF in Eq[2Ql An analytical expression 
of this PDF has been given in (jAime fc Soummerll2004af) . 



The m ean and variance of the intensity can be obtained by several ways. A first method ([Goodman 
1975L I20Q0I ) . is to express the mean intensity E[J] and the second order moment of the intensity E[/ 2 ] as 
a function of C(r) and S(r). The second order moment for the intensity is the fourth order moment for 
the complex amplitude: E[/ 2 ] = E[(C + S)(C* + S*)] 2 ] (omitting the variables r for clarity), which can be 
simplified using the properties of Gaussian distributions, with E[55*55*] = 2E[55*]E[55*] = 2I 2 we obtain: 
E[/ 2 ] = I 2 + 4I C I S + 2I 2 . A second method is to derive a general analytical expression for the mom ents of the 
Ricia n distribution. This can be be obtained either from the definition of t he moments o f Eq |20l ( Goodman] 



19751 ). or computing the derivatives of the moment generating function (jAime fc Soummer 
instantaneous intensity in the focal plane (Eq |15p can be written as: 

I = \C{r)\ 2 + \S{r)\ 2 + 2Re[C*(r) S(r)}. 



2004al ). The 



(22) 



Since E[5(r)*] = E[5(r)]* = (circular Gaussian distribution), the mean intensity is simply the sum of the 
deterministic diffraction pattern with a halo produced by the average of the speckles: I c + I s or I c + / s , 
respectively for direct and coronagraphic images. The variance also finds a simple analytical expression, and 



1 The Mathematica software llWolframlll999h can be used to derive these expressions, and the equivalence between Eq |19l and 
Eq 1201 can be verified easily using the functions Simplify and FunctionExpand. 



0.02 



Fig. 7. — Probability Density Function of the light intensity at 3 different positions in the focal plane, 
corresponding to different amplitude C(r) (or intensity levels I c ). The width of the distribution clearly 
increases with an increase in the level of the constant intensity background. This approach provides an 
alternative explanation of speckle pinning, where the constant background corresponding to the perfect part 
of the wave amplifies speckle fluctuations. 
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we have: 

E[7] = I s +I c 
aj = I 2 S +2I S I C . (23) 

The variance associated with photodetection can be added to this expression to obtain the total variance 
a 2 = a 2 + op, where op is the variance associated to the poisson statistics: op = I c + I s . 
The total variance is therefore: 

o 2 =I 2 + 2I S I C + I C + I S . (24) 



In the case of direct images, the term I c corresponds to the perfect Point Spread Function (PSF) scaled 
to the SR. In the case of coronagraphic images, the focal plane intensity is not invariant by translation, and 
therefore it is technically not a true PSF. However, we will use the term "coronagraphic PSF" for simplicity 
and to follow the general usage in the community. The term I s = E[|5'(r)| 2 ] is a function of the radial distance 
r, which can describe an actual AO halo. Th ese PSFs and halo structures have been studied analytically 



( Moffat 19691 : Racine! 1996 : Racine et al. 1999h . It is also possible to determine the halo profile directly from 
numerical simulations, and an illustration of I c and I s is shown in Fig[8l The long exposure PSF profile is 
the sum of these two contributions, the halo clearing effect for higher SR (jSivaramakrishnan et al.ll200l[) is 
clearly visible between the tw o figures. The shape of the halo is due to the spatially filtered wavefront sensor 
( Povneer fc Macintoshl 2004 ) used in this simulation. In Figl9]we show the effect of a coronagraph on the I c 



-5.5 



-6.5 



-7.5 




-12.5 



-15. 
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Fig. 8. — Numerical simulation to illustrate the decomposition of the mean intensity PSF into two compo- 
nents I c and I s for two Strehl Ratios 90% (V=8) and 95% (V=4) for a direct, non-coronagraphic image. 
The I s term corresponds to the mean speckle halo and the I c term corresponds to the perfect PSF, scaled to 
the Strehl Ratio, so that the total intensity remains normalized (the difference between the two I c profiles 
is neglected here in log scale). The simulation is made with PAOLA. 

term, while the I s term is left unmodified as explained in Sec(21 In this figure we only consider one of the 
previous two AO cases. In this example, the coronagraph is good enough to render the constant background 
term I c negligible when compared to the speckle term I s . 
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Fig. 9. — Effect of a coronagraph on the I c term while the I s term is assumed unmodified (Sec|2|). We only 
consider one case of AO in this figure (SR=95%) and illustrate the I c term for the direct and coronagraphic 
case. The effect of the APLC coronagraph here is to reduce the perfect PSF below the speckle halo. The 
corresponding long exposure image is totally dominated by the halo and no residual ringing remains. 



3.1.3. Effect of a coronagraph on speckle pinning 



In a direct, non coronagraphic image, the term coupling the deterministic C(r) a nd random parts S(r) 



in Eg l22l corresponds to the s o -called "speckle pinning" , discussed by several authors (IBloemhof et al 



Sivaramakrishnan et al 



2002 



Blocmhoi 



2003 



Sivaramakrishnan et al 



200J 



Perrin et al. 



. 2003; 



2001 



Bloemhoi 



2004bh . using an expansion of |.F[P(r)e- , ' ¥, ( I ')]| 2 . In this expansion approach, the PSF consist of a sum of 
terms, some of which contain a multiplicative factor T[P(r)], corresponding to the diffracted field (the Airy 
amplitude for a circular aperture). These terms inherit the zeros of T[P(r)] and contribute all together 
to the pinned speckles. Pinned speckles therefore have zero amplitude at the Airy nulls, an d non-pinned 



speckl es remain at these locations. The first order PSF expansion term, denoted by pi in iPerrin et al 



(J2003J) , corresponds to the pinned speckles and the second order (p2haio) to non-pinned speckles. Higher 
order terms contribute to pinned and non-pinned speckles. This expansion approach provides particularly 
i nteresting insight into the spatial properties and symmetries of the speckles for each order of the expansion 
|Sivaramakrishnan et Zll2003l : [Perrin et alJEiol lBloemholl2004bh . The expansion and our decomposition 
are therefore very similar. In Eq|22j our pinning term 2Re[C*(r) S(r)) includes the contribution of all pinning 
terms from the infinite expansion. However, our constant term in the decomposition is not 1 and this term 
can include the effect of static aberrations, as discussed below. 

The analysis of the statistical properties of the speckles enables a deeper understanding of the pinning 
phenomenon. As shown in Eqj22] and Eqj23j the pinned speckle term of Eq|22] does not contribute to the 
mean intensity, but only contributes to the variance. A numerical simulation is used in FigJT] to illustrate 
the Rician distribution for a direct image, at three different positions in the field: one at the top of an Airy 
ring (strong pinning effect), one at a PSF zero (no speckle pinning) and one at an intermediate position. 
Speckle pinning can be well illustrated by the analysis of these PDFs, as speckle intensity and fluctuations 
are amplified by the term I c . This can be seen in the PDFs in Fig [71 where the widths increase with I c . 
Depending on the amplitude of the Airy pattern at successive rings, the intensity I c is alternatively large and 
small and the variance of the speckles is amplified accordingly by the coherent part of the wave amplitude 
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C(r), with corresponding intensity I c . At the zeroes of the PSF, no amplification occurs and the statistics 
is equivalent to that of a fully developed speckle pattern (exponential statistics). It is important to note 
that speckles fluctuations are not fully cancelled there (FigJTJ), but simply not amplified and their statistics 
is that of laser speckles. 

The effect of a coronagraph on speckles can be well explained from this statistical modeling. Formally, 
as shown in Eqfrjl the effect of a coronagraph is to replace the telescope PSF I c by the on-axis coronagraphic 
PSF after a coronagraph I c . This term includes the effect of static aberrations if they are included in the 
model. 



In the perfect case of an ideal coronagraph achieving a to tal extinction of the star ( Rouan et al. 2000; 
Aime et al.lbooj iKuchner fc Traublbooa Ifoo et al.lEooj) , speckle pinning is fully canceled if there 



arc 



no static aberrations in the system. 



In the case of static aberrations in the pupil (due to polishing and alignment errors for example), 
th ese aberrations propagate thro ugh the coronagraph, according to the description given in SecfSl or 



Sivaramakrishnan et al.l (|20071 ). The deterministic response of the coronagraph I c (Eq[TT]) therefore 



includes the effect of these static aberrations. As a result, static aberrations leaking through a coro- 
nagraph contribute to speckle pinning, even in the case of a perfe ct coronagraph. Such effect ca n be 
produced for example by dead actuators on the deformable mirror (jSivaramakrishnan et al.ll2006l ). 



Following Aime fc Soummer ( 2004bh . we can illustrate the effect of a coronagraph on speckle noise, by 
breaking down the total variance (Eg HM)) into two contributions, a 2 = 2I C I S +I C and a 2 — I 2 +I S . The term 
a 2 is the part of the variance that can be removed by a coronagraph by changing I c , thus affecting speckle 
pinning. In FigfTUland FigfTTlwe illustrate these variance contributions for direct and coronagraphic images. 
In the case of the direct image, the ringing aspect of the variance corresponds to speckle pinning (where the 
variance is amplified) . The coronagraphic variance profile is smooth (no amplification of the speckles is seen 
at the location of the diffraction rings). Details of the numerical simulation are given in Sec l4.2t 



3.1.4- Test of the Rician distribution with numerical simulations 
Tests of the Rician s t atisti cs on real data using the Lick observatory AO system have been carried out 



bv lFitzgerald fc Graham! (|2006f ). showing that the Rician model is consistent with the data. Several compli- 
cations exist in real data, so we test the Rician distribution on numerical simulations to determine whether 
the model is acceptable in a simple case consistent with the hypothesis used in the physical model, with- 
out any additional complicating circumstances or noise. We used PAOLA to generate 10000 AO-corrected 
instantaneous phase screens corresponding to an ExAO system, on a 8m telescope (we used a telescope 
geometry compatible with Gemini or VLT). The parameters chosen for this simulation include 44 actuators 
across the pupil, an integration time and time lag of 0.5ms for a magnitude V = 8 star, and observations in 
the H-band. The atmosphere include a typical C 2 profile for Cerro-Pachon and the seeing is assumed to be 
1.4 arcec. The Strehl Ratio of these simulated images is 83%. 

In each image, we extract the intensities values along a radius to construct 50 intensity series in the 
focal plane at these 50 pixel locations. An example of the first 200 intensity values at an arbitrary location 
is given in FigfT2l 



For each of these 50 points, we performed a Maximum Likelihood estimation of the parameters I c and 
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Fig. 10. — Total variance for direct imaging (blue solid line) with contributions of the variance terms a s (red 
dashed line) and a c which correspond to the pinning contribution (black dotted line) . In direct imaging, the 
variance budget is totally dominated by the speckle pinning effect, at least close to the axis. This simulation 
is made for a SR = 95%. 
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Fig. 11. — Total variance for coronagraphic imaging (blue solid line), with contributions of the variance 
terms a s (red dashed line) and er c which correspond to the pinning contribution (black dotted line). A 
coronagraph can remove the speckle pinning contribution to the variance of the noise. This simulation is 
made for a SR — 95% and using and Apodized Pupil Lyot Coronagraph. 
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Fig. 12. — Example of instantaneous PSF simulated with PAOLA, including both phase and amplitude errors 
(left) with the intensities for 200 independent frames at an arbitrary position in the focal plane (right). 10000 
independent PSFs have been generated at a sampling of 3 times Nyquist, and the values along a radial axis 
are stored for this study. 



I s , assuming a Rician distribution for the data. The Likelihood has been computed as function of the two 
parameters for the un-binned data and maximized using optimization routines of Mathematica. We then 
perform the x 2 an d Kolomogorov-Smirnoff test statistics on these results. We use ten identical identical bins 
for the x 2 test. In Fig[T3J we show two examples of binned data with error bars (here due to the Poisson 
statistics), superimposed with the Rician distribution for the estimated I s and I c parameters. 




Bin number Bin number 



Fig. 13. — Examples of fitted Rician distributions to the binned data, in two cases corresponding to two 
different locations in the focal plane: At a top of an Airy ring (left) and at a zero of the perfect PSF where 
the distribution is exponential (right). 1000 independent realizations have been used in this simulation. 



Our implementation of the Kolmogorov-Smirnoff test is based on a Monte-Carlo estimation of the 
distribution. Indeed, when parameters are estimated from the data (here I c and J s ), the distribution of 
the KS test is not known analytically and must be estimated from Monte-Carlo simulations. For that, we 
calculated a distribution of KS test values for 100,000 random drawings of random Rician data where the 
parameters have to be estimated from the data. This empirical distribution of the KS test values is used to 
generate the empirical right tail values for the test statistic. 

The results of the two tests (right tail values) are given in FigQTj for each of the 50 points along the 
radial axis (between and 8.3A/D). It is interesting to note that for the first two points (on-axis intensity 
distribution), both tests conclude that the Rician model is incompatible with the data at the 5% level. 
Although the statistical properties at the central point are not directly relevant for evaluating the detection 
limits in high contrast images, this questio n is important in itse l f as it corresponds to the distribution of Strehl 
Ratio. This problem has been studied bv lSoummer fc Ferraril (|2007l ) where it is shown that the distribution 
at the center of the ima ge is given by a "rev e rsed" no n -central Gamma distr i bution . This problem was 
studied independently bv lGladvsz et alj (|2006h : Idadyszl (|2006h ; Ichristou et"al] |2006h with similar results. 
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Outside of the center, most pixels pass the two tests except for about 2 pixel locations, which are rejected, 
by one test or the other. This result is consistent with the 5% level we chose, given the total number of 
locations (50 pixels) and we conclude that the Rician model is compatible with the simulated data. We 
have verified on a few sets of simulations that the location of these pixel failing the test is not relevant and 
simply due to the statistics. On the contrary, the central pixel fails the tests systematically. It is interesting 
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Fig. 14. — x 2 an d KS test results for the 50 locations along a radius in the focal plane. The first point 
corresponds to the on-axis case, and the last point is located at 8.3A/-D, with a separation of A/(6Z?) 
between each point. The solid line shows the 5% rejection threshold. The first two points are rejected by 
both tests, which illustrates the non-Rician distribution. In the rest of the field both tests are consistent 
with the Rician distribution hypothesis at the 5% rejection level. 



to compare the estimated parameters with the actual parameters known from the simulation. In Fig |15l 
we compare the actual I c profile known from the simulation data and compare it to the reconstructed I c 
profile from the estimation of the Rician statistics. Both curves show a very good agreement, even somewhat 
surprisingly at the center w here the Rician model i s wro ng. This result is compatible with the satisfactory 
results of the test statistics. Fitzgerald fc Graham ( 2006^ discussed this type of approach as a possible way 
to reconstruct a telescope perfect PSF from the statistics of the PSFs. They note that this procedure might 
be difficult to implement with real data, but is at least verified in simulations in this simple case. 
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Fig. 15. — Comparison between the actual perfect term I c from the simulation without wavefront errors 
(red dashed lines), and the I c values retrieved from the Maximum Likelihood estimation of the parameters 
(solid line), assuming the Rician distribution for the intensity. Both curves show a very good agreement, 
consistent with the test statistics results. Log intensities are plotted as a function of the radial position in 
pixels for easier comparison with other figures. 



3.2. Effect of quasi-static aberrations on the statistical properties 



quasi-static speckles ( 


Beuzit et al. 


1997; 


Oroenheimer et al. 


2001: 


Marois et al. 


2003; 


Boccaletti et al. 


2003, 


2004; 


Soummer et al. 


2006). Practical solutions have been proposed to overcome the quasi-static speckles 



Marois et al.ll2006aj ) o r considered for spac e projects (jBorde fc Traub 



2006 



Give'on et al. 



2006). In ExAO 



coronagraphic images, iHinklev et a"D |2006h have studied residual speckle lifetimes and evidenced two types 
of speckles, one with lifetimes of a few seconds and the other with lifetimes of a few hundred seconds. 
Although there is no reason to believe the universality of these particular values, this problem of quasi-static 
speckles is general enough to need further understanding and theoretical insight. Quasi-static speckles are 
produced by the slow deformations of the optics (thermal and mechanical) that happen for example as the 
telescope slews to follow the star. The typical time scale for these slow variations of aberrations are tens of 
minutes. 

We can generalize the approach of SecfS] to include two types of aberrations with different time scales. 
We will assume that the error wavefront at the entrance pupil consists of the coherent addition of two 
terms: a fast-varying wavefront related to atmospheric residuals, and a quasi-static aberration term. We 
assume that quasi-static aberrations can be described as the sum of a deterministic static aberration and a 
slowly evolving, zero-mean complex aberration: the static aberrations correspond to the polishing errors and 
aberrations produced by the optical design and actual misalignments in the system. Quasi-static aberrations 
corresponds to slow zero-mean fluctuations around this static (and therefore deterministic) contribution. 
Under these hypothesis, we can add a fourth term to the wavefront decomposition of Eqf2] 



#i(r) =A + A s (r) + ai(r)+a 2 (r), 



(25) 



where, as previously, A corresponds to the perfect wave and ^4 s (r) to the static aberrations (deterministic). 
The error terms ai(r) and 02 (r) correspond respectively to the AO-corected and quasi-static aberrations. 
We assume that the lifetime of ai(r) and a 2 (r) are respectively t\ and r 2 > T\. We consider the longest 
lifetime t 2 as the time unit, and we denote by N the number of fast-speckle realizations during a slow-speckle 
lifetime: 

n w 2 

where Ni and iV 2 are the number of realizations for each processes during a long exposure of duration At. 
Following the results of the previous sections, we can write the field amplitude in the focal plane during N 
successive short time intervals T\ as: 



* 4 (r) = <7(r) + Sf{v) + S 2 (r), k = 1, 



(27) 



As previously, C(r) = C(r) + S s (r) where C(r) corresponds to the focal plane field (direct of coronagraphic) 
in the perfect case without aberrations, and S s (r) to the static aberration (e.g. polishing errors or actual 
optical design aberrations), ^'(r) corresponds to the AO-corrected atmospheric aberrations and S^r) to 
the quasi-static aberrations. Outside of the central point in the focal plane, £{ fc ' (r) and S 2 (r) are independent 
circular Gaussian distributions with zero mean: S [ f ] (r) ~ Af c (0,I al ) and S 2 (r) ~ A/" c (0,7 s2 ) with respective 
lifetimes t\ and t 2 . Also, the spatial properties of ^'(r) and 5 2 (r) are different as they result from 
different aberrations sources, and their corresponding variances I a \ and I s i can be generated considering the 
appropriate spatial power spectra. 
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During a time unit t 2 corresponding to the lifetime of a slow speckle, the intensity is: 



/ = ]r|C(r) + sf ] (r) + S 2 (r)| 2 

k=l 
N 

= El^l 2 - ( 28 ) 

k=l 

The expected value of \Uk\ 2 is E[|[4| 2 ] = var[C4] + |E[[/fc]| 2 , where Uf. follows the Gaussian distribution 
Uk ~ J\f c (C(r), I s2 + I a i). Therefore, the mean intensity is: 

E[I}=N(I s2 +I sl +I c ). (29) 

Reminding that I c — \C(r)\ 2 . This expression is consistent with the simple case with one type of speckle 
(Eal23p. as it is expressed here for an exposure t 2 = Nt±. 

The variance of the intensity is defined as: 

N 

var[/] = £var[|tf fe | 2 ]+£cov[|tf fe | 2 ,|E/,| 2 ] 

fc=l k^l 

= Nvir[\U h \ 2 ]+N{N-l)cov[\U k \ a ,\Ui\ 2 ]. (30) 

The covariance cov[|J7fe| 2 , \Ui\ 2 ] can be easily computed using the exp ansion of the high order mome nts of 
a complex Gaussian vector as a function of its high order cumulants (jMcCullag hl ll987t lFerrar l l200fih . The 
covariance finds a simple expression: 

cov[|[/ fc | 2 , \U t \ 2 } - |cov[C4, U t *}\ 2 + E[[7 fe *]E[C/Jcov[[7 fe) C/f] + E[U k ]E[U?]cm[U lt Ui]. (31) 

With Uk ~ Af c (C(r),I S 2 + I s i) and k ^ I cov{Uk, U*] = I S 2, we obtain immediately the covariance and the 
variance (for k=l): 

covOt/fel 2 ,^! 2 ] =I 2 . 2 + 2I c I s2 
var[|C/ fc | 2 ] = {I s2 + hi) 2 + 2I c (I s2 + I sl ) (32) 

Finally, the variance of the intensity fEq|30p. for an exposure t 2 , becomes: 

var[7] = N(I 2 1 +NI 2 . 2 + 2I c (I sl +NI s2 ) + 2I sl I s2 ) 

= Na 2 , (33) 

where the variance cr| for a short exposure T\\ 

g 2 = I 2 sl + NI 2 2 + 2I c (I sl + NI s2 ) + 2I sl I s2 (34) 

is consistent with the case with one type of speckle fEql23|). if I s2 = 0. 

• In this generalized expression of the variance, the pinning term 2I c (I s i + NI s2 ), still exist. However, 
it now includes a term corresponding to the pinning of quasi-static speckle. The coefficient N (the 
ratio of the speckle lifetimes) can be very large, which makes this pinning particularly efficient: a 
quasi-static halo I s2 which is N times lower than I s \ produces the same pinning effect in direct images. 
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This pinning contribution is still directly affected by a coronagraph, which reduces (or cancels) the I c 
term. In the case of a real, imperfect coronagraph, although no atmospheric pinning might be present, 
residual pinning of quasi-static speckles may occur due to the amplification by the factor iV if the 
quasi-static aberrations are not sufficiently small. 

• The cross term 2I s iI S 2 corresponds to speckle pinning of the atmospheric by quasi-static aberrations 
(coherent amplification of the atmospheric speckles by the quasi-static speckles). If the system is 
dominated by the quasi-static speckles (I S 2 large), this term is negligible compared to the halo NI 2 2 . 



It is possible to specify the performance of a coronagraph to avoid speckle pinning, using this expression for 
the variance. Pinning terms are negligible compared to the halo terms if: 



2I c (/ sl +iV/ s2 )«/f 1 +iV/ s 2 



(35) 



assuming the cross term negligible. In the case of only one type of speckles, or when the syst em is dominated 
by quasi-static speckles, this condition simply becomes I c <C I s (|Aime fe; Soummerll20Q4br i. The necessary 
rejection of the coronagraph can be defi ned as a ratio of t he coronagraphic PSF I^ oro to the direct PSF 
jdirect _ gi m ii ar approach was used by iBloemhofl (|2004al ) to quantify the coronagraph effect on pinned 
speckles using the expansion approach. 



4. Dynamic range in coronagraphic images 



In this paper we adopt a definition of the dynamic range based on an expression of the Signal to Noise 
Ratio (S/R), assuming a typical 5a level detection limit. Based on the knowledge of the noise statistics, more 
advanced method s are possible, where the dynamic range can be defined in terms of probability of detection 



and false alarm ( Michel fe; Ferrari 2003 : Ferrari et al. 20061 ). In particular, because of the non-Gaussian 



statistics (exponential distribution for a perfect coronagraph), confidence levels have to be carefully defined, 
as they do not correspond to the usual values for a Gaussian distribution, C. Marois et al. (2007, in 
preparation). In this section we define analytically the dynamic range in a coronagraphic experiment, using 
the variance expressions obtained above. This approach enables a semi-analytical method. We present some 
simulations results for the dynamic range. 



4.1. Signal to Noise ratio and dynamic range 

For simplicity in the calculation of the variance fSec l3.2|l . we expressed the intensity for a time unit r 2 as 
a sum of short exposures t\ (Eq |2"5)) . In the expression of the variance we derived (Eq |3"3")l , the terms I s \, I s2 , 
I c correspond implicitly to a number of photons during an exposure t%. However, it is convenient in practice 
to define normalized intensity terms for a single photon at the entrance aperture and scale them using the 
star flux F+. A normalization term t 2 F* is therefore necessary in order to use such normalized intensities. 
Recalling that a long exposure consists of N2 exposures of durations r 2 , we obtain the final expression of the 
variance for a long exposure At: 



Speckle = N 2 T i F * var[i] 

= JVi t 2 F 2 {I 2 sl + NI 2 s2 + 2I c {I sl + NI s2 ) + 2J sl J s2 ). 



(36) 
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2. 



The photon noise is obtained by applying the same normalization, multiplying Eg 1291 by A^r-j* 

photon = {h + I si + /.a) At, (37) 

and introducing a p = I c + I s i + 7 S 2 for consistency with the case with one type of speckle. The signal 
from a companion is F p f(9) At, where F p is the planet flux and f(8) is a function that describes the off-axis 
response of the coronagraph to a companion. f(9) can be computed independently. For simplicity, we use the 
maximum intensit y of the off-axis comp anion for the signal, but other approaches such as matched filtering 



npE 

can also be used (jSoummer et al.ll2006l ). The signal to noise ratio (S/N) is: 



S/N= , MB^ (38) 
iNnfFtoj + FtAta* 



This expression for the signal to noise ratio can be immediately converted into an expression for the dynamic 
range, assuming a given S/N level. The dynamic range, denoted d is simply the intensity ratio between an 
off-axis companion and the star that produces a given S/N level: F p = dF+. Since At = Nin, we have: 

(39) 




f(0) V A< V 1 nF* ■ 

Other noise contributions can be easily added in this expression if necessary, but a detailed study of any 
specific instrument is beyond the scope of this paper. 



4.2. Semi-analytical method based on the statistical model 

The dynamic range fEq |39[) . can be computed directly from numerical simulations of the terms I s \, I S 2, 
and I c . We follow the equations used in the model (Eq[3]to Eqfl2j) to construct numerical estimates for I c 
and I s (or J a i, I S 2, if applicable). The term A is obtained from its definition (Eq[3|). We consider the average 



over the aperture and over a few independent realizations, generated with PAOLA (jJolissaint et al.ll2006h 
and including both amplitude p(r) and phase screens <p(r). As described in Sec[2j the term \A\ 2 corresponds 
to the Strehl Ratio, and is used to generate the I c profiles: direct and coronagraphic normalized PSFs are 
calculated in the perfect case without aberrations, and multiplied by \A\ 2 . 

The term S(r) is generated following the notation of EqlTOl S(r) = J-[a(r) P'(r)], where a(r) = 
p(r) exp Iv ^ r ) —A is the zero-mean complex amplitude at the aperture. The I s term is obtained averag- 
ing independent realizations of ^(r)! 2 . It is easy to verify that the total intensity I = I c + I s remains 
correctly normalized with this procedure. 

Finally, radial profiles are generated by averaging these results azimuthally. Satisfactory profiles are 
obtained with only a few realizations (typically a few tens), where purely numerical simulations require large 
numbers of independent realizations to estimate the same variance. Once the various I c and I s terms have 
been generated, they can be readily combined according to Eqj36]or Eq[38]to produce variance or detection 
plots. Multiple instrumental cases can be easily studied this way, combining direct or various coronagraphic 
images with different cases of AO is done rapidly using our analytical expressions (for example to study 
the effect of stellar magnitude on the dynamic range). This method also enables the study of quasi-static 
speckles, which is otherwise computationally expensive. 

We compared the variance obtained with the analytical model with the variance obtained from a purely 
numerical simulation, computing a large number (~ 1000) of independent realizations. This process is time- 
consuming because of the wavefront propagation through the coronagraph for each phase screen. Excellent 
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agreement with a full numerical simulation of the variance is obtained with a few realizations (~ 20), using 
the semi-analytical method fFigJTBj) . We note that the model gives an incorrect estimate of the variance 
behind the focal plane occulting mask. This is expected given the approximations described in Eq llll and 
Eq|121 This is not a problem when studying dynamic range, as these expressions are weighted by the off-axis 
transmission f(6) of the coronagraph, and results are only relevant outside the inner working angle of the 
instrument. 
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Fig. 16. — Validation of the results obtained with the semi-analytical model, compared to a purely numerical 
experiment. The variance obtained from numerical estimation of I c and / s based on 20 realizations of phase 
and amplitude screens. Solid black line: model, dashed red line: variance for 1000 realizations, dotted blue 
line: variance estimated spatially on a single frame. 
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4.3. Results and comparison with numerical simulations 



We present some dynamic range simulation results, using the semi- analytical method based on Eql39l 
We consider an eig ht-meter telescope, an AO sy stem with 44 actuators across the aperture, a spatially filtered 
wavefront sensor (jPovneer fc Macintosh! 120041 ) . and an open loop frequency of 2.5 kHz. The atmospheric 
simulation uses typical seeing values and profiles at Mauna Kea. A few cases have been simulated, 
with stellar magnitudes ranging from V=4 to V=8. Strehl Ratios are obtained between 87% and 95%. 
At mospher i c scint illation is also simulated as amplitude screens with PAOLA, which implements the method 



of iRoddierl (|l98ll ). We assumed an instrumental throughput of 25% and an arbitrary exposure time of 
At = 1000s. The atmospheric lifetime is assumed to be 40ms for observations in the H-band and the stellar 
magnitude is V= 4. The coronagraph in the simulation is an APLC similar to the one under study for GPI 



rtSoummerll20"05h . With a 5(7 detection level for this simulation, we illustrate the dynamic range (including 



photon and speckle noise only) in Fig[17] In this example, photon and speckle noise are approximately at 
the same level inside the AO control region. Outside the control region, the dynamic range is set by the 
speckle noise level. 



It is interesting to note that the relative position of the speckle and photon contributions do not depend 
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on exposure time (Eq |39|) . Indeed, only the speckle lifetime and stellar flux have an effect on the relative 
contributions of speckle and photon noise to the dynamic range. This is also true when both fast and quasi- 
static speckles are present. In terms of the contributions of a 2 and a 2 c to the dynamic range as defined in 
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Fig. 17. — Dynamic range simulation for a 1000s exposure (black solid line), showing the contributions of 
the photon noise alone (blue dotted line) and speckle noise alone (red dashed line). In this example, speckle 
and photon noise are balanced inside the control region and the outer region is dominated by the speckle 
noise. 



Aime fc Soummerl (]2004bf ). the APLC coronagraph is almost perfect in this simulation and the part of the 
variance corresponding to speckle pinning, a 2 , has been completely removed. The case illustrated in Fig llOl 
and FigQl] corresponds to the same simulation and shows the effect of the APLC coronagraph canceling the 
pinning contribution completely. 

Finally, we performed a simulation including both fast and quasi-static aberrations. We considered for 
that a simple approach where the power spectrum of the quasi-static aberrations follows a power law /~ 3 , 
with a lifetime of 600s (lOmin). This choice is consistent with the assumption that quasi-static aberrations 
result from polishing errors, usually well described by such power laws. We considered both photon and 
speckle noise, according to Eq |39l In this example, although the quasi-static wavefront error is very small 
(lOnm RMS), it dominates the error budget and limits the dynamic range inside the control radius of the 
AO system. This result depends directly on the power spectrum chosen to generate these aberrations. This 
could explain the results obtained with the Lyot project coronagraph where the dy namic range plots do not 
show the expected halo-clearing region within the control radius of the AO system Elinklev et al. ( 2006h . A 
slightly shallower power spectrum, or a higher level of quasi-static aberrations would easily reproduce the 
observed dynamic range by f illing the cleared region comp letely. This effect could be due mainly to the 
presence of broken actuators (jSivaramakrishnan et al.ll2006l ) in the AO system. 



4.4. Discussion on the effects of speckle reduction techniques 

The noise budgets expressed in Eq l36l and Eq l37l cnable a general understanding of high contrast imaging 
for exoplanet detection. In the case of direct imaging without coronagraph, the main limitation is set by 
speckle pinning which appears as an amplification of the variance levels at the position of the maxima of the 
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Fig. 18. — Dynamic range illustration (black solid line) showing the contribution of fast residual atmospheric 
errors (red dashed line) and quasi-static aberrations (blue dotted line). The simulation includes both photon 
and speckle noise, and the quasi-static aberrations are assumed to follow a f~ 3 power law with a 600s 
lifetime. 



diffracted pattern. In the absence of static aberrations, pinning happens at the top of the Airy rings. In the 
presence of static aberrations and even with a perfect coronagraph, residual pinning exists at the position 
of the maxima of the coronagraphic PSF corresponding to the propagated static aberrations through the 
coronagraph. For example, in the case of low-order static aberrations described by Zernike polynomials, 
the location of the maxima of the PSF (direct or coronagraphic) dictates where pinning happens. These 
locations are no longer at the top of the perfect Airy pattern in this case. Quasi-static aberrations also add 
to the problem by creating additional pinning terms fEq |34|) . This phenomenon was observed at the AEOS 
telescope with the Lyot pr o ject coronagraph, where pinne d speckles have been a ssociated with dead a ctuators 
on the deformable mirror ( Sivaramakrishnan et al. 2006h . or AO waffle mode ( Makidon et al. 2005h . 



In reality, coronagraphs are not perfect and suffer from various sources of errors, including chromaticity, 
manufacturing imperfections, alignment errors, etc. The cancellation of speckle pinning in the final error 
budget provides a constraint on the coronagraph performance. This can be achieved based on the comparison 
between the terms I c and / s described above. 

With current state of the art AO systems and coronagraphs alone, it is clear that the residual noise 
terms (even if the pinning term is well controlled) would not enable the detection of giant planets with typical 
contrasts of 10~ 7 at 0.5 arcsec, or Earth-li ke planets with contrast of 10~ 10 at 0.1 ar csec. Additional speckle 
reduction techniques have been developed ( Borde fc TraublEood : Give'on et al ]l2006h to improve the contrast 
performance further in combination with wavefront control and coronagraphs. These techniques affect the 
term I s by modifying the pupil aberrations to create a dark zone in the search region of the fie ld. In the case 



of ground-based high contrast imagers, for example with the GPI speckle calibration system (jWallace et al 



2006), the speckle reduction will only affect the I S 2 term, as the sensing involved in the technique is slow. 
Some speckle reduction techniques may also affect the residual speckle pinning in some region of the field, 
as th ey create a static wavefront aberration w hich creates a particularly asymmetric PSF presenting a dark 
zone (jCodona et alj|2006t ISerabvn et al.ll2006l ). 



The combination of coronagraphy and speckle nulling will therefore remove most of the terms contribut- 
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ing to the noise variance aj . The limit for the association of coronagraphy and speckle nulling is set by the 
level of residual fast atmospheric aberrations. Also, by decreasing the term J S 2, the speckle reduction system 
puts a stronger constraint on the coronagraph, whose performance must suppress pinned speckles to a level 
below other speckles in the dark region produced by speckle nulling. 

Finally, the remaining parts of the speckle noise which have not been removed by a coronagraph or by 
the speckle reduction system, can be removed further using post-processing differential methods, combinin g 
multi- waveleng t h info r mation from dual-ba nd imaging, or an integral field spectrograph Marois et al.l (|2000) ; 



Sparks fc Fordl (|2002f k iMarois et al.l (|2005l ). or using differential rotation iMarois et al 



(|2006ah . These dif- 



ferential methods affect all the residual speckle variance contributions (pinning and halo) in the same way. 
The combination of these three stages (coronagraph, speckle reduction and speckle processing) is necessary 
to reach the ne cessary dynamic range regime for exoplanet imaging. This approach has b een chosen for the 
design of GPI (jMacintosh et al.l l2006h . and have been implemented in the laboratory bv iTrauger fc Traubl 
d2007h . 



5. Conclusion 

This paper develops theoretical insight in high contrast imaging, using a statistical model to understand 
the performance limitation s of coronagraphic experiments. We discuss the application of the statistical model 
( Aime fe Soummer 2004bh to the case of coronagraphic images, and show that the statistical model is valid 



outside of the focal plane mask occulting area. We show some limitations of the statistical mo d el at the 



central point of the image, and a modification of the model is proposed by ISoummer fc Ferraril (|2007l ) to 
solve this problem. 

It has been confirmed by several observations that the main limitation for high contrast imaging is 
due to the presence of quasi-static speckles. These do not average out over time, and are difficult to 
calibrate. We presented the first theoretical attempt to understand the effect of quasi-static speckles and their 
interaction with other aberrations (static and residual atmospheric). We obtained an analytical expression 
for the variance of the intensity, which includes speckle and photon noise in the presence of static, quasi- 
static and rapidly-varying aberrations. This result enables the use of a semi-analytical method, which has 
been compared successfully with purely numerical simulations. This semi-analytical method enables fast 
simulations and the possibility to compare and combine easily various types of parameters for instrument 
design studies. It also enables a breakdown of the noise into specific contributions (pinning, halo, speckle, 
photon). This method can be used to model real observations and extract information on relative noise 
contributions from various sources (residual atmospheric noise, quasi-static noise). 

The model also provides an understanding of the speckle pinning effect in high Strehl images, where 
bright speckles appear at the location of the diffraction rings. This effect, explained by a coherent ampli- 
fication of the speckle noise, and can be cancelled by a coronagraph. The theory provides insight on the 
required level of performance for the coronagraph to effectively cancel speckle pinning. In the presence of 
quasi-static aberrations, we showed that pinning also exists between the perfect part of the wave and quasi- 
static aberrations, but is weighted by the ratio of the lifetimes (which can be a large number of the order of 
a hundred or more). This puts strong requirements on the level of static and quasi-static aberrations and 
may explain the results currently obtained with the Lyot project coronagraph, the first instrument of its 
kind to use an extreme AO system. Detail modeling of the Lyot project data, and specific simulations for 
the future generation of high contrast imager (GPI), will be done in the future. 
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This model was developed for the monochromatic case. Generalization to wide band is not straight- 
forward, as the speckle noise is highly correlated with wavelength, which makes this generalization difficult, 
but can be used for multi- wavelength speckle subtractions. Also, a generalization of the model to the case of 
subtraction residuals in multi-wavelength imaging would be particularly interesting and is a topic for future 
study. 

Although this presentation is focused on the ground-based case, the same formalism can be directly 
used for space-based observations, with similar results, but different orders of magnitudes (in particular for 
the speckle lifetimes) . An identical expression of the variance of the noise was established independently by 



Shaklan et all (|2005l) for the TPF error budget. 



Once a coronagraph has removed the speckle pinning contribution, the dynamic range can be improved 
further by use of speckle reduction techniques (speckle nulling) or speckle subtraction (ADI, multi-wavelength 
or polarization) . The effect of speckle nulling on the noise budget is now well understood from the expression 
of the variance. The model can be used to draw the requirements on the relative performance of the 
coronagraph and speckle reduction system. This theory brings an understanding of the effects on the dynamic 
range of the various elements in a high contrast instrument, including wavefront control, coronagraphy, 
speckle reduction and differential calibration techniques. All these stages are crucial and necessary to reach 
high contrast imaging, and they are going to be implemented in the next generation of high contrast imagers, 
currently in development. 
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